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mallest.  Force  vs.  vertical  displacement  r.nd  force  vs.  horizontal  strain  at 
he  center  of  the  disc  were  plotted  for  each  test  specimen.  No  appreciable 
volume  effect  on  the  tensile  strength  of  the  specimens  was  apparent. 

A  one-dimensional  progressive  failure  model  is  proposed,  consisting  of 
inearly  elastic  springs  in  parallel  which  can  deform  independently  of  each 
>ther.  The  springs  are  assumed  to  be  perfectly  brittle  and  the  values  of 
stiffness  and  strength  are  assigned  using  the  Monte  Carlo  technique.  While 
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mmerous  variables  involved. 

A  two-dimensional  progressive  failure  finite  element  model  is  presented.) 
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PREFACE 


This  report  covers  the  first  year  accomplishments  in  the 
research  program  entitled,  "Mechanical  Behavior  of  Rock  Under 
Static  Loading,  "  R.  W.  Heins  -Co- Principal  Investigator.  The 
program  is  Part  A  of  the  project  entitled,  "Aspects  of  Mechanical 
Behavior  of  Rock  Under  Static  and  Cyclic  Loading"  (Contract 
H0210004).  Part  B  of  the  project  is  published  in  a  separate 
volume. 
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ASPECTS  OF  MECHANICAL  BEHAVIOR  OF  ROCK  UNDER  STATIC  LOADING 

PART  A 


Summary  of  Work  to  Date 

Valders  limestone  was  tested  for  indirect  tensile  strength  using  the 
well  known  Brazilian  test  Various  sizes  of  samples  were  tested.  The 
largest  sample  was  three  inches  in  diameter  by  two  inches  long  (3"  x  2"), 
and  the  smallest  was  one  inch  in  diameter  by  one-half  inch  long  (1"  x  1/2"). 
Results  of  these  tests  are  summarized  in  Table  1.  4.  1.  Force  versus  vertical 
displacement  and  force  versus  horizontal  strain  at  the  center  of  the  disc 
were  plotted  for  each  test  specimen.  A  closed-loop  servocontrolled  machine 
(MB  machine)  was  used  for  the  test.  Application  of  load  was  controlled  by 
lateral  strain  rate  at  the  center  of  the  disc.  Complete  stress-strain  (force- 
displacement)  curves  were  obtained. 

A  one-dimensional  progressive  failure  model  is  given  in  Chapter  2. 

This  model  consists  of  linearly  elastic  springs  in  parallel.  The  springs  de¬ 
form  independently  of  each  other.  They  are  assumed  to  be  perfectly  brittle, 
i.  e. ,  any  additional  load  (force  or  displacement  loading)  beyond  the  peak 
will  cause  a  sudden  failure  of  the  spring.  The  spring  stiffness  and  strength 
(failure  displacement)  values  are  assigned  to  the  elements  by  the  Monte  Carlo 
technique,  using  various  types  of  probability  distribution  functions.  Five 
types  of  distribution  functions,  Uniform,  Left-skewed,  Right-skewed,  Uni¬ 
model  and  Bimodel,  are  used.  The  one-dimensional  model  computer  program 
is  given  in  Appendix  A.  The  following  variables  were  analyzed:  (1)  number 
of  springs,  (2)  strength  variance,  (3)  strength  distribution,  (4)  stiffness 
variance,  and  (5)  stiffness  distribution. 

A  two-dimensional  progressive  failure  finite  element  model  is  given 
in  Chapter  3.  In  this  model,  the  Brazilian  test  was  modeled  by  a 


two-dimensional  finite  element  structure.  Elastic  properties  and  strength 
values  for  each  element  of  the  finite  element  model  are  assigned  by  Monte 
Carlo  technique  using  a  uniform  distribution  function.  A  series  of  elastic 
analyses  were  performed  on  this  structure.  In  each  cycle,  the  element  with 
worst  combination  of  principal  stress  and  strength  is  found.  Ninety-nine 
percent  of  this  element  stiffness  is  then  substracted  from  the  entire  structure 
stiffness.  The  elastic  analysis  is  performed  on  this  new  structure,  and  the 
cycles  are  repeated  until  catastrophic  failure  occurs.  The  following  variables 
are  analyzed:  (1)  ratio  of  allowable  compressive  strength  to  tensile  strength, 

(2)  number  of  loaded  node  points,  (3)  strength  variance,  and  (4)  stiffness 
variance. 

Conclusions 

The  following  conclusions  may  be  extrapolated  from  the  limited  amount 
of  work  presented  in  this  report. 

(1)  In  order  to  obtain  a  complete  force-displacement  curve  for  brittle 
materials,  such  as  Valders  limestone,  the  loading  rate  must  be  controlled  by 
horizontal  displacement  (or  strain)  rate  in  order  to  avoid  catastrophic  failure 
at  or  very  near  the  peak  load. 

(2)  The  variation  in  the  volume  of  the  Brazilian  test  samples  did  not 
appear  to  have  any  significant  effect  upon  the  tensile  strength  of  Valders 
limestone.  The  variation  in  average  tensile  strength  was  3.  7%  around  the 
mean  value. 

(3)  The  one-dimensional  model  does  not  have  any  direct  use  in  prac¬ 
tice,  but  it  serves  a  very  important  use  in  this  study.  It  serves  as  a  simple 
model  to  study  numerous  variables  involved. 

(4)  The  strength  variance  and  the  shape  of  the  strength  distribution 
has  a  marked  effect  upon  the  force-displacement  curves.  On  the  other  hand, 
the  stiffness  variance  and  the  shape  of  the  stiffness  distribution  does  not 
have  a  very  significant  effect  upon  the  force-displacement  curves. 
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(5)  The  use  of  bounded  distribution  functions;  allows  the  effect  of 
variables  to  be  examined  with  physical  interpretations  to  be  placed  on  these 
bounds. 

(6)  While  the  two-dimensional  model  gives  good  qualitative  results, 
it  should  not  be  considered  a  quantitative  tool  for  solving  rock  mechanics 
problems  until  further  modifications  and  improvisations  are  made. 

(7)  The  two-dimensional  model  is  capable  of  showing  how  failure 
develops  in  mixed  stress  fields  such  as  Brazilian  tests. 

Future  Work 

The  general  direction  of  the  future  program  will  be  as  follows: 

(1)  Additional  Brazilian  tests  on  Valders  limestone  will  be  performed. 

(2)  Two  other  rock  types,  Dacite  and  St.  Cloud  Gray  granodiortte, 
will  be  also  used  for  Brazilian  tests. 

(3)  Additional  work  on  the  two-dimensional  model  will  be  performed. 
Work  is  underway  on  a  two-dimensional  model  using  isoparametric  finite 
elements,  automatic  mesh  generation,  and  other  refinements.  Anisotropy  of 
elements  will  also  be  included. 

(4)  Other  types  of  two-dimensional  models  are  also  planned.  The 
greatest  emphasis  will  be  placed  upon  new  failure  criterion.  Literature 
search  is  underway. 

(5)  After  a  new,  improvised,  two-dimensional  model  is  developed, 
work  on  a  three-dimensional  model  will  be  started. 
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CHAPTER  1 


BRAZILIAN  TEST  ON  VALDERS  LIMESTONE 


1.  1  Introduction 

Valders  limestone  was  tesiod  for  indirect  tensile  strength  using  the 
well-known  Brazilian  test  Various  sizes  of  samples  were  tested.  The 
largest  sample  was  three  inches  in  diameter  by  two  inches  long  (3"  x  2"), 
and  the  smallest  was  one  inch  in  diameter  by  one-half  inch  long  (1"  x  1/2"). 
Force  versus  vertical  displacement  and  force  versus  horizontal  strain  at 
the  center  of  the  disc  were  plotted  for  each  test  specimen.  A  closed-loop 
servocontrolled  machine  (MB  machine)  was  used  for  the  test.  Application 
of  the  load  was  controlled  by  lateral  strain  rate  at  the  center  of  the  disc. 

The  primary  objective  of  the  laboratory  testing  program  was  to  obtain 
complete  stress-strain  (or  force-displacement)  curves  for  qualitative  compar¬ 
ison  with  the  mathematical  models  proposed  in  Chapters  2  and  3.  The  ob¬ 
jective  behind  using  various  sizes  of  samples  was  to  study  the  effect  of 
the  sample  volume  on  the  tensile  strength  as  determined  by  the  Brazilian 
test  It  was  also  intended  the .  the  mathematical  model  would  give  similar 
qualitative  information  regarding  the  effect  of  the  sample  volume. 

1.  2  Brazilian  Test 

The  Brazilian  test  for  measuring  the  tensile  strength  of  brittle  mater¬ 
ials,  like  concrete,  rocks,  etc. ,  is  not  only  widely  used  in  this  country  but 
in  Europe  and  many  other  parts  of  the  world  also.  It  is  probably  the  simplest 
test  in  terms  of  preparing  the  samples,  carrying  out  the  testing,  and  analyz¬ 
ing  the  results.  A  circular  disc  is  cut  from  the  material  to  be  tested  and  is 
loaded  by  a  theoretically  concentrated  line  load,  Figure  1.  2. 1.  The  stress 


o 
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FIGURE  1.  2. 1  Brazilian  Test 


along  the  load  axis  in  x  direction  is  constant  up  to  about  +  0.  8  R  from  the 
center,  and  is  given  by  Timoshenko  and  Goodier(l)*  to  be 


c  -  ZA  .  2L. 

t  7rR  TTDt 


(1.  2.  1) 


where 


St  =  tensile  stress  perpendicular  to  the  load  axis 
4  =  concentrated  line  load 

R  =  radius  of  the  disc 
t  =  thickness  of  the  disc,  and 
D  =  diameter  of  the  disc. 


The  dimensions  of  the  specimen  can  be  measured  before  the  testing 
and  the  load  P  needed  to  fail  the  sample  is  recorded  from  the  loading  machine. 
Then,  Equation  1.  2.  1  gives  the  tensile  strength  of  the  specimen. 

The  tensile  strength  of  the  specimen  is  obtained  by  Equation  1.  2.  1  if 
the  failure  is  caused  by  excess  tension  along  the  load  axis.  However,  many 
times  the  failure  is  initiated  by  a  "crushing  and  wedging  mechanism"  at  the 
point  of  application  of  the  load  P.  Figure  1.  2. 1.  Coates  (2)  gives  an  empir¬ 
ical  relationship  in  order  to  take  this  mechanism  into  account  Fairhurst  (3) 

* 

Number  in  parenthesis  indicates  reference  number  in  the  bibliography  at 
the  end  of  this  chapter. 
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has  suggested  that  this  "crushing  and  wedging"  action  can  be  avoided  by 
applying  the  load  over  an  arc  of  about  7°. 

In  the  present  laboratory  investigation  considering  the  samples 
to  date,  the  crushing  was  confined  to  a  very  small  portion  of  the  loading 
surface,  thus  Equation  1.  2.  1  was  usod  to  calculate  the  tensile  strength* 

1.  3  Sample  Size 

The  sample  sizes  were  decided  based  upon  the  following  criteria: 

1.  The  volume  of  the  sample  should  be  doubled  each  time,  for  the 
same  diameter  sample. 

2.  The  failure  load  of  the  sample  should  not  exceed  the  capacity  of 
the  machine  (25  kips). 

3.  The  diameter  of  the  sample  should  be  such  that  it  fits  into  the 
vertical  displacement  measurement  device,  Figure  1.  4.  2, 
conveniently  and  easily. 

The  third  criterion  limited  the  samples  to  be  three  inches  in  diameter 
or  less.  The  second  criterion  limited  the  largest  sample  size  to  be  three 
Inches  in  diameter  by  two  Inches  in  thickness  (3"  x  2").  Therefore,  begin¬ 
ning  with  l  sample  3"  x  2",  the  following  sample  sizes  were  selected  to 
meet  the  first  criterion. 

1.  3"  diameter  x  2",  1"  and  1/2"  thick, 

2.  2"  diameter  x  2",  1"  and  1/2"  thick,  and 

3.  1"  diameter  x  2'',  1"  and  1/2"  thick. 

1.  4  Test  Procedure  and  Results 

A  closed-loop  sei*vocon trolled  machine  (MB  machine),  Figure  1.  4.  1, 
was  used  for  the  Brazilian  test  of  Valders  limestone.  Application  of  load  was 
controlled  by  the  horizontal  strain  rate  at  the  center  of  the  sample  disc. 
Applied  force  versus  horizontal  strain  at  the  center  output  was  directly  re¬ 
corded  by  the  MB  machine  recorder, and  a  special  additional  plotter  was  used' 
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Figure  1.  4. 1 


The  MB  Universal  Testing  Machine 
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(foreground,  Figure  1.  4. 1)  to  plot  force  versus  vertical  displacement. 

A  special  platen  was  constructed  to  measure  the  vertical  displace¬ 
ment,  Figure  1.  4.  2.  The  top  and  bottom  plates  were  very  carefully  machined 
to  make  them  perfectly  horizontal  such  that  no  eccentricity  will  be  Introduced 
in  the  application  of  the  vertical  load.  An  adapter,  to  fit  the  platen  at  top 
Into  the  load,  was  also  very  precisely  made  in  order  to  have  perfect,  full 
contact  with  the  load  cell.  One  strain  gage  was  mounted  on  top  of  each  of 
the  two  cantilevers  extending  from  the  top  platen.  First,  the  cantllevers  are 
puiled  down  by  tightening  the  nuts  on  top  of  the  two  poles.  As  the  bottom 
platen  moves  up,  the  cantilever  will  be  displaced  upward  and  will  record  a 
certain  value  of  strain.  The  platen  set  was  calibrated  to  convert  strain  on 
the  cantilever  to  vertical  movement  of  the  platen. 

To  measure  horizontal  displacement  on  the  axis  perpendicular  to  the 
applied  load,  another  displacement  measuring  mechanism  was  devised, 

Figure  1.  4.  3.  A  considerable  amount  of  effort  was  spent  on  this  device  to 
perfect  it  in  an  attempt  to  make  it  give  satisfactory  results;  however,  satis¬ 
factory  results  were  not  obtainable  and  the  device  was  discarded. 

Because  of  sudden  failure  and  permanent  deformation  of  the  device, 
it  was  decided  to  use  strain  gages  directly  mounted  at  the  center  of  the  rock, 
Figure  1.  4.  4.  Two  SR-4  strain  gages  were  mounted  on  the  sample,  one  each 
on  the  two  faces  of  the  disc,  Figure  1.  4.  4-  This  alternative  is  expensive 
because  two  strain  gages  are  expended  for  each  sample  tested,  and  is  time 
consuming  because  considerable  time  Is  spent  in  mounting  the  strain  gages. 
However,  the  strain  gage  method  gives  an  accurate  strain  value  at  the  center 
of  the  disc.  The  average  of  the  two  strain  gage  readings  was  used  in  the 
force  versus  horizontal  strain  plots. 

After  the  samples  were  cored  and  cut  to  proper  sizes,  they  were  dried 
in  an  oven  for  48  hours.  Two  strain  gages,  one  on  each  side,  were  then 
mounted  in  the  standard  manner.  The  sample  was  positioned  in  the  loading 
frame  and  carefully  centered,  Figure  1.  4.  5. 
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Figure  1. 4.  2  Vertical  Displacement  Measurement  Device 


Figure  1.  4.  3  Horizontal  Displacement  Measurement  Device 
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Figure  1.  4.  4  Brazilian  Test  Sample 
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Figure  1.4.5  Brazilian  Test 
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The  loading  rate  in  the  first  phase  of  the  testing  program  was  con¬ 
trolled  by  the  vertical  displacement.  But  we  encountered  the  same  problem 
as  Hudson,  et  al.  (4)  did.  The  rock  samples  broke  very  violently  (shattered 
to  pieces)  at  the  peak  load.  This,  obviously,  did  not  satisfy  our  objective 
of  obtaining  a  complete  stress-strain  carve,  and  another  method  of  obtaining 
a  complete  stress-strain  curve  must  be  found.  Hudson,  et  al.  were  success¬ 
ful  in  using  horizontal  strain  rate  at  the  center  of  the  disc  to  control  the 
loading  on  the  sample  and  to  obtain  a  complete  stress-strain  curve.  There¬ 
fore,  in  this  laboratory  investigation  we  also  decided  to  use  horizontal 
strain  rate  to  control  the  load  on  the  rock  sample. 

The  MB  machine,  Figure  1.  4. 1,  was  programmed  to  reach  the  highest 
value  of  strain  allowed  by  the  machine  in  800  seconds.  The  amount  of  ver¬ 
tical  loading  exerted  by  the  machine  was  programmed  to  be  proportional  to 
the  horizontal  strain  rate.  If  the  horizontal  strain  rate  exceeded  the  preset 
(programmed)  value  in  the  machine,  then  the  vertical  load  will  decrease, 
and  vice-versa.  With  this  servocontrolled  machine,  it  is  thought  that  the 
change  in  vertical  loading,  depending  upon  the  horizontal  strain  rate, 
occurs  within  a  few  microseconds. 

There  were,  on  the  average,  three  samples  of  each  size  tested,  see 
Table  1.  4.  1.  A  typical  complete  vertical  force-vertical  displacement  plot 
is  shown  in  Figure  1.  4.  6,  and  a  vertical  force-horizontal  strain  plot  is  shown 
in  Figure  1.  4.  7. 

The  testing  stopped  either  when  the  MB  machine  went  into  oscillation* 
and/or  when  the  strain  gages  mounted  on  the  sample  broke.  As  the  vertical 
crack  along  the  loaded  axis  opened  up,  it  also  cracked  the  strain  gages 
securely  mounted  on  the  sample,  Figures  1.  4.  9,  1.  4. 10  and  1.  4.  11.  Approx¬ 
imately  half  of  the  tests  stopped  due  to  failure  of  the  strain  gage.  The  other 
half  of  the  tests  were  terminated  when  the  machine  went  into  oscillations, 

P 

It  is  believed  that  the  machine  went  into  oscillation  when  the  horizontal 
strain  rate  was  fluctuating,  forcing  the  closed-loop  servocontrolled  machine 
to  load  end/or  unload  too  rapidly. 
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Vertical  Di  spla cement  1"  =  590x  10 


Horizontal  Strain  at  Center 
1"  =  200  x  10~6  inch/inch 
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Figure  1.4.10  Sample  7M  After  Failure 


19 


see  Table  1.  4.  1.  In  general,  the  samples  with  smaller  volume  failed  due  to 
strain  gage  failure,  Figures  1.  4.  9,  1.  4.  10  and  1.  4. 11,  whereas  for  larger 
volume  samples  testing  was  stopped  when  the  machine  went  into  oscillation, 
Figure  1.  4.  8.  Figures  1.  4.  8,  1.  4.  9,  1.  4.  10  and  1.  4.  11  also  show  the  var-* 
latlon  from  test  to  test  in  the  type  and  amount  of  fracture.  As  was  expected  (1), 
due  to  highly  localized  compressive  stress  at  the  point  of  application  of  load, 
the  failure  crack  in  general  propagated  from  the  outside  to  the  center  of  the 
disc  along  the  loaded  axis. 

1.  5  Discussion  of  Results 

Vertical  applied  force  versus  vertical  displacement  and  vertical  applied 
force  versus  horizontal  strain  at  the  cente  of  the  disc  plcti,  like  Figures 
1.  4.  6  and  1.  4.  7,  were  obtained  for  al.  the  samples  tested.  Qualitatively, 
the  plots  for  all  the  samples  were  similar. 

After  the  initial  almost  linear  increase  in  force  versus  the  displacement 
for  both  types  of  plots,  the  peak  was  achieved  in  about  50  to  80  seconds. 

The  time  to  reach  maximum  horizontal  strain  was  set  at  800  seconds.  For 
about  three  to  four  seconds  after  the  peak  was  reached,  there  was  no  notice¬ 
able  movement  of  the  plotter  pens.  At  this  time,  crushing  at  the  load  points 
and  the  propagation  of  a  hairline  crack  to  the  center  of  the  disc  was  evident. 

The  rock  disc  did  not  fail  suddenly  at  this  point  because  the  loading  was 
controlled  by  horizontal  strain  rate  and  not  by  the  vertical  force  or  displace¬ 
ment  rate.  If  the  loading  was  controlled  by  the  vertical  force  or  displacement 
rate,  then  it  is  obvious  from  Figure  1.  4.  6  that  the  sample  would  have  failed 
suddenly  at  the  peak  load — catastrophic  failure.  The  use  of  the  horizontal 
strain  rate  to  control  the  machine  prevented  the  catastrophic  failure. 

Immediately  after  the  peak  point  was  reached,  the  crack  weakened 
the  sample  and  reduced  the  load  bearing  capacity  and  the  unloading  occurred. 

The  sample  was  still  carrying  the  reduced  load,  and  ths  horizontal  strain  in 
the  weakened  sample  increased  (the  crack  widened).  This  process  continued 
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until  the  test  was  stopped  due  to  breakage  of  the  strain  gage, or  due  to 
machine  oscillations, or  when  the  maximum  strain  value  was  reached  in  the 
preset  total  time  of  800  seconds. 

Table  1.  4. 1  gives  some  details  regarding  the  Brazilian  tests.  The 
mean  of  the  average  tensile  strength  was  1360  psi,  Variation  in  the  average 
tensile  strength  thus  was  +  50  psi  or  3.  7%  around  the  mean  value.  From 
these  limited  tests,  we  can  then  conclude  that  the  sample  volume  did  not 
have  much  effect  on  the  tensile  strength  of  the  Valders  limestone. 

Variation  in  tensile  strength  within  tests  was  of  the  higher  magnitude, 
Table  1.  4.  1.  For  all  the  one-inch  diameter  samples,  the  variation  was  16.  3%, 
for  two-inch  diameter  samples  it  was  6.  2%,  and  for  three-inch  diameter 
samples  it  was  9.  7%.  The  mean  of  the  within-tests  variations,  for  all  the 
samples,  was  10.  7%.  Thus  the  within-tests  variation  in  tensile  strength  is 
significantly  higher  when  compared  to  the  variation  in  tensile  strength  due 
to  volume  change. 

1. 6  Conclusions 

From  the  limited  amount  of  laboratory  testing  completed  so  far,  the 
following  conclusions  may  be  extrapolated: 

1.  In  order  to  obtain  a  complete  force-displacement  curve  for  brittle 
materials,  such  as  Valders  limestone,  the  loading  rate  must  be 
controlled  by  horizontal  displacement  (or  strain)  rate  in  order  to 
avoid  catastrophic  failure  at  or  very  near  the  peak  load. 

2.  The  variation  in  the  volume  of  the  samples  did  not  appear  to  have 
any  significant  effect  upon  the  tensile  strength  of  the  Valders  lime¬ 
stone.  The  variation  in  average  tensile  strength  wa&  3.  7%  around 
the  mean  value. 
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1.  7  Future  Program 

The  general  direction  of  the  future  program  will  be  as  follows: 

1.  More  samples  of  Valders  limestone  of  each  size  will  be  tested,  using 
stress  (load)  controlled  loading  for  comparative  purposes. 

2.  Two  other  rock  types  will  be  tested,  with  both  strain  rate  and  stress- 
controlled  loading. 

3.  Statistical  methods  will  be  used  to  aid  in  evaluation  of  these  results. 
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CHAPTER  2 

ONE-DIMENSIONAL  MODEL 


2.  1  Introduction 

A  one-dimensional  progressive  failure  model  is  proposed  in  this 
chapter.  This  model  consists  of  linearly  elastic  springs  in  parallel.  The 
springs  deform  independently  of  each  other.  They  are  assumed  to  be  per¬ 
fectly  brittle,  i.  e. ,  any  additional  load  (force  or  displacement  loading) 
beyond  the  peak  will  cause  a  sudden  failure  of  the  spring.  The  spring  stiff¬ 
ness  and  strength  (failure  displacement)  values  are  assigned  to  the  elements 
by  the  Monte  Carlo  technique,  using  various  types  of  probability  distribution 
functions. 

t 

2.  2  Distribution  Functions 

Five  type*  of  distribution  functions  are  used  to  randomly  generate 
the  spring  stiffness  and  strength.  They  are;  unifor  n,  left-skewed,  right- 
skewed,  unimodal  and  bimodal.  The  uniform  distribution  is  well-known, 
and  the  others  are  derived  from  it  as  suggested  by  Hemmerle  (1).*  These 
distribution  functions  are  chosen  because  they  represent  wide  variations  and 
are  inexpensive  to  generate  on  a  digital  computer.  There  is  no  reason  to 
expect  that  these  distributions  represent  the  properties  of  a  real  rock.  It  is 
hypothesized  that  the  inhomogeneous  real  rock  properties  are  statistically 
inhomogeneous.  The  varieties  of  distribution  functions  should  be  sufficient 
to  examine  their  effect  upon  the  progressive  failure  of  the  assumed  statisti¬ 
cally  inhomogeneous  rock.  An  added  benefit  is  that  the  distributions  are  all 
bounded  over  a  range,  a  to  b,  thus  allowing  physical  interpretation  to  these 
bounds.  The  width  of  each  distribution  in  terms  of  variance,  V,  and  mean 

* 

Number  in  the  parenthesis  indicates  reference  number  in  the  bibliography  at 
the  end  of  this  chapter. 
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value  of  the  distribution  is  tabulated  in  Table  2.  2. 1. 

TABLE  2.  2.  1 


!  Wo. 

Distribution 

Subroutine  Name* 

Mean  Value 

VO 

Width  of  the 
Distributions  =  b-a 

1 

Uniform 

UNIFRM 

a  «■  b 

2 

(12  V) 1/2 

2 

Left-Skewed 

RSKEW** 

2b  +  a 

3 

(18  V) 1/2 

3 

Right-Skewed 

LSKEW 

b  +  2a 

3 

(18V) 1/2 

4 

Unimcdal 

UNIMOD 

a  +  b 

2 

(24V) 1/2 

5 

Bimodal 

31  MOD 

a  +  b 

2 

(8V)1/2 

- - 1 

** 


See  Appendix  A. 

* 

Subroutine  RSKEW  actually  produces  a  distribution  which  is  left-skewed,  i.  e. , 
the  mean  of  the  distribution  is  to  the  right  of  the  midpoint  with  the  tail  of  the 
distribution  stretched  out  to  the  left.  Similarly,  LSKEW  produces  a  right- 
skewed  distribution. 


We  can  see  from  Table  2.  2.  1  that  by  specifying  the  mean  value  and 
variance,  V,  we  can  easily  calculate  the  range  of  the  distribution  (values  of 
a  and  b). 

2.  3  Progressive  Failure  Modal 

The  parallel  spring  model  is  loaded  with  displacement  loading.  The 
total  stiffness  of  the  model  is  first  computed.  It  is  equal  to  the  sum  of  stiff¬ 
nesses  of  all  the  springs.  After  applying  the  first  displacement,  we  can  cal¬ 
culate  the  force  necessary  to  cause  the  displacement  It  is  equal  to  the 
total  stiffness  of  the  model  multiplied  by  the  displacement  This  gives  us 
one  point  for  the  force  versus  displacement  plot.  Before  increasing  the  dis¬ 
placement,  a  check  is  made  as  to  whether  or  not  strength  (specified  in  terms 
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of  displacement)  of  any  spring  is  exceeded.  The  spring  whose  strength  is 
equal  to  the  applied  load  (both  in  terms  of  displacement)  will  fail  at  the  load 
and  is  removed  from  the  model.  The  new  total  stiffness  of  the  model  is 
computed  by  subtracting  the  failed  spring  stiffness.  The  displacement  is 
increased,  the  force  needed  to  cause  the  displacement  is  calculated,  the 
failed  spring  is  picked  and  eliminated,  and  the  cycle  is  repeated  until  all 
springs  have  progressively  failed.  All  the  force  versus  displacement  points 
are  plotted  to  obtain  a  force  versus  displacement  plot  for  the  model. 

2.  4  Discussion  of  Results 

To  reduce  the  number  of  independent  variables  in  the  model  and  facil¬ 
itate  comparison,  a  hypothetical  rock  was  selected  for  modeling  in  a  uniaxial 
compression  test.  The  specimen  is  a  cylinder  with  a  diameter  of  two  inches 
and  a  length  of  four  inches.  The  rock  has  a  modulus  of  elasticity  of  7.  6  x 
106  psi,  giving  tot A  stiffness,  K  =  —-,  of  6  x  106  lb/inch. 

The  one-dimensional  model  computer  program  is  given  in  Appendix  A. 

It  was  used  to  obtain  force  versus  displacement  plots  for  the  model.  There 
are  various  independent  variables  which  need  to  be  analyzed.  It  was  decided 
to  analyze  the  following  variables: 

(1)  Number  of  springs 

(2)  (a)  Strength  variance  (varied,  holding  yield  pcint  constant), 

(b)  Strength  variance  (constant,  yield  point  varied). 

(3)  Strength  distribution 

(4)  Stiffness  variance 

(5)  Stiffness  distribution 

The  effect  of  varying  number  of  springs,  N,  is  of  two  types.  The 
least  important  but  most  obvious,  Figure  2.  4.  1,  is  that  the  plotting  of  smaller 
numbers  of  points  by  a  "connect  the  dots"  rule  will  not  produce  a  smooth 
curve.  More  important,  however,  is  the  statistical  "law  of  large  numbers," 
which  says  that  if  "sufficient"  number  of  samples  is  taken,  then  the  sample 


25 


distribution  is  arbitrarily  close  to  the  theoretical  distributioa  This  means 
that  as  N  approaches  infinity,  the  curve  produced  by  the  model  becomes 
realistic.  Thus  the  number  of  elements  should  be  selected  large  enough  to 
make  the  force-displacement  curve  smooth  and  reproducible,  even  with  larger 
number  of  elements.  However,  increasing  N  will  increase  the  cost  of  analysis, 
in  terms  of  increased  computer  time.  Figures  2.  4. 1,  2.  4.  2,  2.  4.  3,  2.  4.  4 

and  2.  4.  5  show  the  result  of  varying  N(N  =  10,  100,  500  and  1000)  with  uni¬ 
form  strength  and  stiffness  distributions.  The  curves  for  N  s  500  and  N  =  1000 
are  smooth  and  nearly  identical,  Figure  2.  4.  5,  while  curves  for  N  =  10  and 
N  =  100  are  rough  and  not  reproducible.  It  was  concluded  that  500  spring 
elements  were  sufficient,  and  this  number  was  used  for  subsequent  analysis 
of  other  variables. 

The  effect  of  strength  variance  on  the  force-displacement  plot  may  be 
studied  in  two  different  ways:  (1)  by  holding  the  yield  point  cpnstant  and 
varying  the  variance,  and  (2)  by  holding  the  variance  constant  and  varying 
the  yield  point 

Let  As  represent  the  lower  bound  of  the  strength  distribution.  Thus  the 

first  element  will  fail  when  the  applied  displacement  equals  A  .  The  model 

s 

is  elastic  until  the  first  element  fails.  Thus,  A  represents  the  yield  point  of 

s 

the  model.  B  represents  the  upper  bound  of  the  strength  distribution.  When 

S 

the  displacement  becomes  B  ,  the  last  element  fails.  Thus  B  represents  the 

s  s 

ultimate  failure  point  of  the  model. 

The  result  of  chenging  variance  of  strength  distribution,  holding 
A  (=  0.  004  inch)  constant,  is  shown  in  Figure  2.  4.  6.  With  A  =  0.  004  inch, 
the  model  simulates  the  hypothetical  rock  which  has  the  yield  point  of  about 
24,  000  lbs.  or  7,  600  psi.  The  strength  distribution  variance  ranged  from 
zero  to  5.  6  x  10  The  stiffness  distribution  was  held  constant  and  both 
distributions  were  uniform.  With  bounded  distributions  (of  the  types  used 
here),  the  role  of  the  strength  variance  is  to  fix  the  upper  and  lower  bounds 
of  the  strength  distribution  if  the  mean  of  the  distribution  is  specified.  Here, 
the  strength  variance  was  varied  holding  the  lower  bound  constant  at  0.  004 
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Figure  2.4.1  E:fec:L  of  N,  Number  of  Spring  Elements 
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Figure  2.4.2  Effect  of  N 
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Figure  2.4.3  Effect  of  N 
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N  =  30,  100,  500,  1000 
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inch.  Thus,  tho  strength  distribution  fixed  the  upper  bound  or  the  naxlir  ji?i 
displacement  at  which  the  specimen  will  fail.  If  Ag  and  B  are  coincident, 
variance  =  0,  then  the  yield,  peak  and  maximum  displacement  occur  at  the 
same  displacement  and  weakest  link  failure  occurs.  This  implies  that  the 
hypothetical  rock  simulated  has  a  homogeneous  strength. 

The  second  way  to  study  the  effect  of  rock  strength  variance  is  to 
vary  the  yield  point  of  the  rock,  but  holding  the  variance  constant.  The 
strength  and  stiffness  variance  was  held  constant  and  the  uniform  distribu¬ 
tion  was  used.  Ag  varied  from  0  to  0.  008  inch.  The  resulting  force-displace¬ 
ment  plots  are  shown  in  Figure  2.  4.  7. 

From  Figures  2.  4.  6  and  2.  4.  7,  we  can  readily  see  that  the  effect  of 
the  strength  variance  on  the  force-displacement  plot  is  significant. 

Figures  2.  4.  8  to  2.  4.  12  show  the  effect  of  the  five  strength  distribu¬ 
tions  on  the  force-displacement  plots.  The  strength  distributions  are  uniform, 
unimodal,  left-skewed,  right-skewed  and  bimodal,  respectively.  The  uni¬ 
form  stiffness  distribution  was  held  constant.  The  strength  and  stiffness 
variance  and  yield  point  were  also  held  constant.  We  can  see  that  the 

strength  distributions  have  marked  effect  upon  the  force-displacement  response 
of  the  model. 

Effect  of  stiffness  variance  is  shown  in  Figure  2.  4.  13.'  The  strength 
variance  and  the  yield  point,  Ag,  were  held  constant.  The  stiffness  variance 
ranged  from  zero  to  ten  with  the  mean  stiffness  fixed  at  6  x  10®/N,  where  N 
is  the  total  number  of  spring  elements.  The  uniform  distribution  was  used 
for  the  strength  and  the  stiffness  variance.  We  can  see  from  Figure  2.  4. 13 
that  the  stiffness  variance  does  not  have  a  very  significant  effect  upon  the 
force-displacement  response  of  the  model. 

Figure  2.  4.  14  shows  the  effect  of  the  five  stiffness  distributions  on 
the  force-displacement  plots.  The  stiffness  distributions  are  uniform,  uni¬ 
modal,  left-skewed,  right-skewed  and  bimodal,  respectively.  The  uniform 
strength  distribution  was  held  constant.  The  strength  and  the  stiffness 


34 


LlU.1  i-UJl  LiXLU JJOJliJJJJLiJ.U.1  .Lt.1 J  UllU-LU-IXiUN 


O 

O 


o 

oo 


o 


o 

Ni- 


O 

CM 


<>) 


in 

Cvl  ' — v 


o 

cm  O 
O 


h 

1 0 

>  <  *’ — 

ilJ 

□ 

O 

a: 

'"CL 

CO 

I — ^ 

o 


m 


(81  GOOD  BGHBd 


Figure  2.4.8  Strength  Distribution  Uniform 
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Figure  2,4.9  Strength  Distribution  Unimodal 


Figure  2 .4 . 10  Strength  Distribution  Left-Skewed 


Figure  2.4.11  Strength  Distribution  Right-Skewed 
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Figure  2.4>12  Strength  Distribution  Bimodal 


41 


variance  and  the  yield  point  were  also  held  constant  We  can  observe  from 
Figure  2.  4.  14  that  the  stiffness  distributions  do  not  have  a  very  significant 
effect  upon  the  force-displacement  response  of  the  model. 

2. 5  Observations 

The  following  observations  may  be  made  regarding  the  one-dimensional 

model. 

(1)  The  model  does  not  have  any  direct  use  in  practice,  but  it  serves  a 
very  important  use  in  this  study.  It  serves  as  a  simple  model  to 
study  numerous  variables  involved.  This  will  help  a  great  deal  in 
the  development  of  two  and  three-dimensional  model. 

(2)  This  simple  one-dimensional  model  is  relatively  very  inexpensive  to 
use  in  terms  of  computer  time  for  analysis  of  variables.  We  can 
expect  a  similar  effect  of  variables  in  case  of  two-  and  three- 
dimensional  model,  thus  eliminating  a  need  for  detailed  analysis  of 
these  variables  in  subsequent  model  developments. 

(3)  The  use  of  bounded  distributions  allow  the  effect  of  variables  to  be 
examined  with  physical  interpretations  to  be  placed  on  these  bounds. 

(4)  Five  hundred  spring  elements  are  sufficient  for  this  model  to  obtain 
smooth  force-displacement  plots. 

(5)  The  minimum  value  of  strength  distribution  defines  the  yield  point  for 
the  model  and  the  maximum  value  defines  the  ulitmate  failure  dis¬ 
placement 

(6)  The  shape  of  the  strength  distribution  has  a  marked  effect  upon  the 
force-displacement  plot  of  the  model. 

(7)  The  stiffness  variance  and  the  shape  of  the  stiffness  distribution  does 
not  have  a  very  significant  effect  upon  the  force-displacement  plot 

of  the  model. 
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CHAPTER  3 

TWO-DIMENSIONAL  MODEL 


3.  1  Introduction 

While  many  papers  exist  on  the  subject  of  inhomogeneity  in  elastic 
and  plastic  problems,  they  treat  mainly  that  class  of  problems  which  we 
have  called  deterministically  inhomogeneous.  What  we  have  called  statis¬ 
tical  inhomogeneity— which  is  a  more  realistic  model  for  polycrystalline 
materials  such  as  rock— cannot  be  treated  by  exact  methods.  Neither  can 
they  be  easily  handled  by  finite  difference  methods,  due  to  the  difficulty 
with  boundary  conditions  between  regions  with  different  properties. 

The  emergence  over  the  last  decade  of  the  finite  element  method  has 
no  less  than  revolutionized  all  of  continuum  mechanics.  Problems  that  were 
previously  considered  Intractable  or  at  least  very  difficult,  are  now  solved 
routinely.  The  method  has  also  seen  much  use  in  rock  mechanics  as  wit¬ 
nessed  by  Increased  volume  of  papers. 

Inhomogeneous  problems  of  all  kinds  have  been  especially  simplified 
by  the  finite  element  method.  Each  element  can  be  assigned  its  own  unique 
properties  with  no  more  difficulty  than  assigning  the  same  property  to  all 
elements. 

This  approach  was  used  by  Su,  Wang,  and  Stefanko  (1)*  to  study  the 
stress  distribution  around  a  mine  opening  (circular  tunnel).  They  assumed 
that  both  Young's  modulus  and  Poisson's  ratio  were  normally  distributed. 

Not  surprisingly,  in  light  of  the  "law  of  large  numbers,  "  they  found  that  the 
mean  of  the  stress  distribution  for  any  region  approached  the  value  it  would 
have  if  the  material  was  assumed  homogeneous  but  that  large  differences  in 
local  stresses  would  occur  when  compared  to  the  material  being  normally 
distributed  elastic  properties. 

Number  in  the  parenthesis  indicates  reference  number  in  the  bibliography  at 
the  end  of  this  chapter. 
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All  other  papers  encountered  In  a  literature  search  referred  to  deter¬ 
ministic  inhomogeneity  which  is,  of  course,  easily  treated  by  the  method. 

It  will  be  assumed  that  the  reader  is  familiar  with  the  finite  element 
method.  Many  excellent  references  and  textbooks  are  available. 

3*  2  The  Model 

A  finite  element  progressive  failure  model  was  developed  and  pro¬ 
grammed  in  FORTRAN  for  solution  by  a  digital  computer.  While  some  portions 
of  existing  programs  from  several  sources  (2)  were  used,  it  was  found  that 
no  existing  program  could  be  satisfactorily  modified.  The  majority  of  work 
under  the  present  contract  has,  therefore,  been  devoted  to  writing  the  pro¬ 
gram  given  in  Appendix  B. 

Basically,  the  model  consists  of  a  series  of  elastic  analyses  on  a 
two-dimensional  finite  element  structure.  In  each  cycle  the  element  is 
found  which  has  the  worst  combination  of  principal  stress  (computed  from  the 
elastic  analysis)  and  strength  (assigned  as  an  element  property  in  the  program). 
This  element  is  then  considered  "failed"  by  the  program.  This  maximum 
principal  stress  criterion  of  failure  is  clearly  an  oversimplification  of  ele¬ 
mental  failure  mechanism  for  rock.  Since  such  failure  mechanisms  are  the 
subject  of  heated  debate  among  experts  in  fracture  mechanics,  the  principal 

stress  criterion  may  be  considered  an  initial  step  toward  an  eventual  more 
realistic  mod&l. 

From,  the  ratio  of  the  maximum  or  minimum  principal  stress  in  the 
critical  element  to  its  compressive  or  tensile  strength,  the  loads  and  dis¬ 
placements  are  scaled  so  that  the  element  is  exactly  failed,  i.  e. ,  the  ratio 
becomes  one.  The  program  then  computes  the  force  needed  to  cause  the  cal¬ 
culated  displacement  and  reports  the  number  of  the  failed  elements,  the 
external  load  causing  the  failure,  and  the  displacement  under  the  load.  Any 

other  force  or  displacement  quantities  needed  can  be  printed  out  at  the  same 
time. 
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Most  of  the  element  st'ffness  is  removed  by  subtracting  99%  of  the 
element  stiffness  of  the  failed  element.  Not  subtracting  the  entire  stiffn  jss 
avoids  division  by  zero  in  the  matrix  solution  process.  The  idea  of  as?  jming 
the  failed  element  to  lose  almost  all  of  its  stiffness  again  is  an  oversimplifi¬ 
cation  and  should  again  be  considered  but  an  initial  step. 

The  initial  force  and  displacement  conditions  are  reset  and  the  cycle 

is  repeated  as  an  elastic  analysis  of  the  remaining  (unfailed)  portion  of  the 
structure. 

It  should  be  noted  that  the  program  is  not  iterative  in  the  usual  sense, 
i.  e. ,  it  does  not  increment  or  monotonically  increase  the  load  (neither  force 
or  displacement  load).  The  elastic  analysis  is  iterated.  In  each  cycle,  a 
new  structure  is  analyzed.  The  new  structure  essentially  omits  (except  for 

1%)  the  element  failed  in  the  previous  cycle.  The  cycles  are  continued  until 
complete  failure  occurs. 

3.  3  Test  of  the  Modfll 

For  Initial  testing  of  the  finite  element  progressive  failure  model,  the 
diametral  compression  or  Brazilian  test  was  chosen.  It  suits  this  purpose 
very  well  since  it  includes  both  tensile  and  compressive  stress  regions. 

The  Brazilian  test  is  accepted  to  be  an  indirect  tensile  strength  test 
for  rock,  and  a  rather  large  proportion  of  reported  tensile  strengths  for  rocks 
are  obtained  by  this  method.  Recent  investigations  (3)  have  shown  that 
failure  rarely,  if  ever,  initiates  at  the  tensile  stressed  region  in  the  center 
of  the  loaded  diameter.  This  was  also  confirmed  by  the  laboratory  investiga¬ 
tion  outlined  in  Chapter  1.  Failure  is  initiated  by  crushing  under  the  loads 
at  the  enas  of  the  diameter  and  immediately  proceeds  inward,  creating  the 
vertical  crack  commonly  seen.  The  compressive  stress  at  the  ends  of  the 
diameter  is  theoretically  infinite  under  a  line  load  but  decreases  as  the  arc 
of  contact  is  increased.  Fairhurst  (4)  showed  that  tensile  failure  is  theoreti¬ 
cally  possible  if  the  compressive  strength/tensile  strength  ratio  was  greater 
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than  about  12  for  an  arc  of  contact  o=  tan"1  1/12,  but  that  for  a  =  tan'1  1/6 
tensile  failure  is  possible  with  a  ratio  of  8.  (Note:  tan_1l/12  is  about 
4.  8°  and  tan"  1/6  is  about  9.  6°. )  Hudson,  et  al.  (3)  found  that  even  for 

arcs  of  contact  exceeding  10°,  no  tensile  failure  of  the  assumed  kind  could 
be  found. 

A  finite  element  mesh  to  model  the  Brazilian  test  is  shown  in 
Figure  3.  3.  1.  It  consists  of  81  nodes  comprising  128  triangular  elements. 

A  structure  with  five  rings  of  elements  (121  nodes,  200  elements)  was  also 
tried  but  could  not  be  run  due  to  computer  size  limitations. 

Constant  strain  triangles  are  also  probably  preferable  to  more  sophis¬ 
ticated  elements.  At  least  they  will  be  easier  to  work  with.  The  program 
size  is  limited  by  die  computer  co.*e  size  and  determined  primarily  by  the 
size  of  the  stiffness  matrix.  For  a  two-dimensional  structure,  this  size  is 
twice  the  total  number  of  nodes.  Constant  strain  triangles  can  be  shown  to 
give  the  largest  number  of  elements  for  a  given  degree  of  accuracy  of  results 
and  thus  for  a  given  program  size  they  provide  the  maximum  information 
about  tne  progression  of  failure.  The  node  numbering  scheme  shown  gives 
a  bandwidth  of  36,  which  is  felt  to  be  a  minimum  for  the  four  ring  structure. 

The  results  of  six  computer  runs  are  presented  in  the  next  section. 

For  five  of  these  runs,  the  specimen  was  loaded  by  a  vertical  downward  dis¬ 
placement  on  node  number  five  while  the  opposite  node  (number  78)  was 
held  fixed  in  both  directions.  This  simulates  not  a  live  load  but  an  arc  of 
contact  of  11-1/4°  (one-half  element  width  in  each  direction).  Because  the 
high  stress  gradients  under  the  load  cannot  be  well  modeled  by  the  constant 
strain  triangles,  this  probably  corresponds  to  a  somewhat  smaller  arc  on 

real  rock.  One  test  (A-3)  had  three  nodes  loaded  and  three  fixed  simulating 
an  arc  of  contact  of  33-3/4°. 

The  first  three  runs  were  with  homogeneous  element  properties.  The 
fourth  and  fifth  runs  were  with  element  strengths  distributed  two  different 
ways,  and  the  sixth  run  had  only  element  Young's  modulus  distributed. 
Poisson's  ratio  was  held  constant  at  0.  25,  but  would  be  varied  in  later  testing. 
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Properties  were  assigned  to  each  element  by  subroutine  CARDS  using  the 
Monte  Carlo  scheme  developed  for  the  one-dimensional  model.  Only  the 
uniform  distribution  was  used  but  the  other  distributions  may  also  be  tried. 

It  is  doubtful  that  they  will  make  any  large  difference  in  the  results. 

Compressive  strength  for  any  element  is  a  constant  multiple  (provided 
as  data)  of  the  stochastically  generated  tensile  strength.  Since  any  element 
can  fail  only  once,  this  does  not  disturb  the  randomness  of  the  strengths. 

3.  4  Test  Results 

The  conditions  of  the  six  tests  are  shown  in  Table  3.  4.  1.  The  Young's 
modulus  was  a  constant  7,  600,  000  psi  except  in  test  A-6  where  its  mean  was 
7.  6  x  10  psi.  The  mean  (or  constant)  tensile  strength  was  chosen  to  give 
a  mean  compressive  strength  of  13,  000  psi. 

Let  us  discuss  the  results  for  each  test  individually. 

In  test  A- 1,  the  specimen  failed  immediately  by  crushing  the  elements 
adjacent  to  the  loaded  node.  With  these  elements  failed  on  both  ends  of  the 
loaded  diameter,  there  was  little  resistance  to  large  displacement  and  the 
test  was  terminated.  This  test  points  up  a  particular  inadequacy  of  the 
model  in  that  all  load  resistance  would  not  be  removed  by  crushout  under  the 
loads  in  a  real  specimen.  The  crushed  zone  would  retain  considerable  stiff¬ 
ness  and  a  crack  would  propagate  from  the  outside  to  the  center  of  the  speci¬ 
men.  The  failure  initiation  of  the  model  appears  entirely  correct,  however. 

With  a  higher  ratio  of  compressive  to  tensile  strength  in  test  A-2, 
failure  initiated  at  the  center  of  the  disc  and  ran  outward.  The  progression 
of  the  failure  is  shown  in  Figure  3.  4.  1.  Ultimate  failure  was  by  crashing 
under  the  .loads  as  before,  but  initiation  and  propagation  were  as  expected 
for  the  higher  ratio  of  tensile  to  compressive  strength.  Force-displacement 
under  the  load  for  test  A-2  is  shown  in  Figure  3.  4.  2.  It  is  easily  seen  that 
either  force  or  displacement  controlled  loading  will  produce  complete  failure 
simultaneously  with  failure  initiation.  After  n  =  11,  the  specimen  has  failed 
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and  further  results  are  meaningless. 

TABLE  3.  4.  1 


BRAZILIAN  TEST  CONDTTT HMR 

Test 

No. 

Ratio 
-  °C/gT 

No.  of 

Loaded  Nodes 

1 

■04 

Inhomogeneity 

A-l 

8 

1 

1625 

Homogeneous 

A- 2 

24 

1 

541.  67 

Homogeneous 

A-3 

8 

3 

1625 

Homogeneous 

A-4 

24 

1 

541.  67 

' 

Tensile  strength  range  from 

41.  67  psi  to  1041.  67  psi 

A- 5 

24 

1 

541.67 

Tensile  strength  range  from  441,  67 
psi  to  641.  67  psi 

A-6 

24 

I 

541.  67 

Young's  modulus  range  from 

1,  600,  000  psi  to  13,  600,  000  psi 

The  three-point  loading  in  test  A- 3  produced  a  tensile  failure  at  a 
compressive  strength-tensile  strength  ratio  of  8.  This  was  to  be  expected 
oscciuse  a  three-point  loading  will  not  create  high  concentration  of  compres¬ 
sive  stress  at  the  point  of  application  of  load.  In  Figure  3.  4.  3,  the  progressive 
pattern  of  failure  is  shown.  Failure  initiated  at  the  center  in  tension  and  pro¬ 
gressed  simultaneously  in  both  directions.  Ultimate  failure  wos  due  to 
wedging  out  of  an  area  underneath  the  load  which  remained  unfailed.  The 
force-displacement  curve  is  shown  in  Figure  3.  4.  4.  Here,  force  loading 
would  produce  instant  failure,  but  displacement  loading  will  jump  from  n  =  1 
to  n  =  15  and  then  to  ultimate  failure. 

Test  A-4  proved  to  be  too  inhomogeneous  and  failure  was  very  scat¬ 
tered  before  ultimate  failure  due  to  crushing  of  the  bottom  load  point  This 
rather  unrealistic  failure  pattern,  shown  in  Figure  3.  4.  5.  perhaps  shows  that 
degree  of  inhomogeneity  might  be  needed  to  be  defined  to  model  a  real  rock. 

The  force-displacement  behavior  in  Figure  3.  4.  6  shows  another  aspect  of 
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Inhomogeneity.  Failure  in  either  force  or  displacement  loading  would  not 
occur  until  17  elements  had  failed. 

Test  A-5  was  similar  to  test  A-4  but  with  greatly  reduced  degree  of 
inhomogeneity  of  the  strength.  In  this  case,  the  failure  proceeded  from  the 
center  outward  but  slightly  to  the  left  of  the  vertical  diameter  as  shown  in 
Figure  3.  4.  7.  The  peak  force  and  displacement  again  occurred  on  the  first 
element  which  fails  indicating  that  failure  will  be  instantaneous  in  either 
force  or  displacement  loading,  Figure  3.  4.  8. 

The  test  of  inhomogeneous  modulus,  test  A-6,  provided  little  new 
information.  Even  though  the  stiffness  was  highly  inhomogeneous,  the  stress 
distribution  appears  Lo  be  similar  to  the  homogeneous  case.  Apparently  in 
the  Brazilian  test  at  least,  the  strength  is  a  more  sensitive  variable.  This 
was  also  shown  to  be  true  in  case  of  the  one-dimensional  model. 

The  curves  of  horizontal  strain  vs.  vertical  force  are  highly  confus¬ 
ing  as  they  loop  back  upon  themselves  as  the  specimen  unloads  with  element 
failure.  A  typical  example  is  shown  in  Figure  3.  4.  11. 

3. 5  Conclusions 

1.  While  the  model  in  its  present  oversimplified  form  gives  good 
qualitative  results,  it  should  not  be  considered  a  quantitative  tool  for  solv¬ 
ing  rock  mechanics  problems  until  further  modifications  and  improvisations 
are  made. 

2.  Strength  appears  to  be  a  more  sensitive  variable  than  stiffness. 
Further  work  is  needed  to  determine  if  this  is  a  function  of  the  model  and  its 
assumptions  or  of  the  Brazilian  test  geometry. 

3.  The  model  is  capable  of  showing  how  failure  develops  in  mixed 
stress  fields.  This  may  ultimately  prove  to  be  its  most  important  application. 
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appendix  a 

ONE-DIMENSIONAL  PROGRAM 

T;>e  following  pages  contain  a  listing  of  the  main  program  and  seven 
subroutines,  which  were  used  to  produce  the  foree-dlsplacement  plots  In 
Chapter  2.  The  programs  contain  reference  to  numerous  subroutines  provided 
by  tne  Madison  Academic  Computing  Center  (MACC)  as  system  library  rou¬ 
tines.  Most  of  these  are  Involved  with  the  graphics  package  (2)*  which 
produ.~ed  the  plots.  Random  numbers  needed  for  Monte  Carlo  technique  were 
obtained  by  MACC  library  routine  RANUN  (3).*  The  operation  of  the  program 
is  described  briefly  below.  While  the  program  as  presented  can  be  run  only 
on  the  Univac  1108  at  MACC.  It  Is  hoped  that  some  Insights  Into  the  model 
can  be  gained  from  their  presentation  here.  It  Is  further  hoped  that  the  dis¬ 
tribution  subroutines  will  be  useful  to  others  Interested  In  simulation. 

The  main  program  first  fills  arrays  containing  element  stiffness  (STIFF) 
and  strength  (DISM)  with  random  variables  from  the  appropriate  distribution 
by  calls  to  GENS  and  GENK.  It  computes  Initial  total  stiffness  (TSTF)  and 
uses  a  MACC  utility  routine  (URSORT)  to  sort  both  arrays  In  ascending  order 
of  strength.  In  displacement  loading,  the  force  (FORCE)  at  each  load  incre¬ 
ment  Is  computed  as  the  product  ol  current  total  stiffness  and  the  current 
displacement  (DISP1).  The  total  stiffness  Is  then  reduced  by  the  stiffness 

of  the  failed  element  The  arrays  (FORCE  and  DISPL)  are  then  scaled  and 
plotted. 

Subroutines  GENK  and  GENS  fill  arrays  STIFF  and  DISM,  respectively, 
with  random  variables  by  calls  to  appropriate  distribution  subroutines.  These 
a  rray-i  tiling  routines  also  compute  sample  statistics  for  comparison  with 
distribution  statistics. 

The  distribution  subroutines  (UNIFRM,  UNIMOD,  RSkSW.  LSKEW, 
BIMOD)  obtain  random  variates  from  the  five  distributer  ■  as  described  in 

Bibliography  at  the  end  of  Chapter  2. 
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Chapter  2. 

A  list  of  input  variables,  $ome  of  which  are  3et  in  the  program  rather 
than  read  from  data  cards,  is  given  in  Table  A-l, 


1.  N 

2.  VK 

3.  VS 

4.  IK 

(1) 

(2) 

(3) 

(4) 

(5) 

5.  IS 

6.  IPLOT 

(1) 

(2) 

(3) 


7.  AS 

8.  NSER 

9.  EK 


TABLE  A-l 
INPUT  VARIABLES 
Total  number  of  spring  elements 
Stiffness  variance 
Strength  variance 

Code  number  for  the  stiffness  distribution 
Uniform  distribution 
Unlmodal  distribution 
Left- skewed  distribution 
Right-skewed  distribution,  and 
Bimodal  distribution. 

Code  number  for  the  strength  distribution.  Same  code  as  IK. 

Code  number  for  force-displacement  plots 
Plot  not  needed 

Separate  plots  needed  (e.  g. ,  Figure  2.  4.  4) 

Superimpose  plots  (e.  g. ,  Figure  2.  4.  5) 

Minimum  strength,  yield  point 

Serial  number  of  the  data  set 

Mean  of  the  stiffness  distribution  equals  6000000.  /N, 
for  ail  distributions. 
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r 


^  3!  N  OR  l!C  PAM,  ,(jNr  r)TM.  #  ...... . 

DTHfNSI'lV  STIFF  n  OO0».niSM(  1000  1 
1  Xt’A  TA  f  1  01  .  Y"A  TA 1 1  0)  ,  XI  2  )  .  v< 

F ou ivsL f«rr  mrsMi  n  .or srl  i2) ) 

CALL  l NIT  PL (1C. 10 .ft) 

Y<1  ):0. 

XI? )T 10. 

V(1  )  =  0. 

Y(? >r loo. 


RICDL(  1000)  .PliPCEl  10C0  >  » 


l  RTAO  101  N.VK  .vs.  T*  .IS  ,  IPLUT.45.MSfi? 
lF(,N.ir.o.)  rn  ru  ioc 

PRINT  10-> 

c  form  ••,rirr\r^'-  matrix  . . 

call  gfnk  (N.TK.VK.STTFF) 
r  F»1»M  F  T  Pf  MG  T  M  V  A  T  5  T  x  . 

pflCL  1  >  r  M  c.  (M.IS.V^.OISM.AS)* . 

C  CUMPUTf  TriTAL  S  T  T  OF  NF  FS  .. 

t<tf*o,  . . 


C 

c 


ru  ?  I  r  1  ,  ,v 

?  TSTF:  TStr  ♦c.tiffiT) 
r3UPT  TN  A  SO  F.\  OINC  URO'-p  u  F  STPFNCTNF 

nrcnCflLL  <JP<f"*T  U.N.  PTSM.STTFF.nuMl.O^OUM?) 

n ifpi. ArtMLM  leaping . 

FRJNT  107  N<;rq 
CIS  °L  11)10. 

FOR  CE  ( 1  i:  0. 
rn  i  i  =  ’.m 

FtlRCE  II  )r  TFTf  *Df«;pL  (|  , 

Tn>TF“T^Tr-5.  TIFF  (  T  —  J  J 
CONTINUE 


Reproduced  from 
best  available  copy. 


f 

4 


our  hr 1,0 SP  CH  1  0,N*  rtwcr  »L  .pfak.c.  .Pin 
P^INT  1  On  L  .F  ORCriL). HI  SPU  l  ).DT<vPL(N 
on  r  i:i.m 


FflRCF  1  ?)  .PTSPL  (?  ) 


HTSPL  II  ):  ni«;PL(  I)  *1000. 

FOR  Of  I]  )r  CtiPCFI  J  )  eQ.ooi 

C.0  TO  I  1.  <i  .  S  )  .IPLflT 

CALL  OPLMT  IKl.NRL0K.r?,m.N4) 

CALL  RA Gr  1 1 1 . . ft. S. O.^HrtO»DFP ) 

CAL  L  SC  AL  F  I  X  .2.X0A  TA  ,  J  ,c.,  5> 

CALL  SCALF  (  Y  .2  .YOATA  ,  1  ,  10.  S) 

CALL  'iASI  r  I  2  •  •  2  »  •  X  P  A  T  A  ,  o.  .O..YOATA.F 
CALL  PLtlT  12.  ,2.. ‘’HUP) 

XMA  Jr  ?. 


00.) 


00  10 

no  IF  J.-.1.4 


XMJ  NrXMAj*|j*o,p7) 

CALL  PLOT  1XMTN, 2.0, 4HCUWM) 
CALL  ’’L  VC  TP  I  o  .  #  ,  i  .  (iHPulilN  ) 

1*  CALL  PLOT  1XMIN.2..2HUP) 

Y«*A  J*  YM  a. 1*1  .  * 

CALL  PLOT  (  XM  A  J  .  7.  « 40'UJWN  ) 
CA||.  Pi  VT  TP  |  O  .  #  0 .  1  r  .  4  mPdua)  ) 
TMUMr s. i 


10 


YNlJMr  XM  .AJ-0 . 1 

CAL  L  PL  MU  MR  |  y  mum.)  .  ?  ,  IN  UM  .  ?H  T  7 . 0  .  I  4 . 0  .  > 
C«LL  'HOT  1  YM  AJM.  f  ?Hup  j 
CALL  PLOT  I  ?  ML!  P  ) 

YMA  j:  ■>.  n 


Of)  20  Jrl.r 


HO  2C  jrl  ,  0 
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Y *1  n:  v*  AJ  ♦  ( J»  n.  i  > 

call  ( r.  . ymt n. imrsiwM 

r«u  r>i  vr  tp  i  o.  1 .  o, .  4ihju«n  > 

->'  CALL  0L:!’  i  .YH’M.  >HJJP  ) 

VMA  j-  YMAJ.  1  . 

CALL  1) T  I  ’.  .YKA  JtUHnuW*  ) 

CALL  ’LV^TR  «  0.1c  .0.  .4H!)UWN> 

T  MU M: ?0*T 
YN'UM:  yiY  AJ  -0.0  7 
TALL  PL  MW  MP  M  .R.  YN1JM,  INUM.JWM.O.  14 ,0.  ) 

■>n  TALI.  PLOT  (  7  •  0 »  YM  A  J  »  2:<UP  1 

CALL  UL  Tr  XT  (  1  .  7 . * u . V URC F  ItOOO  I.  *  > t §  •  *  0.  ?1  •  *0  .  » 

C*LL  PL  Tr  X  T  <  Y.7C  .  1 .4»  •  mSPLACFMRNT  (.001  IN' .  I  »  .  0.  ?  1  , 0  .  ) 

C  At  L  l  !Nr  (^1  *PL»  XT  AT  A  . FU7CF . YDA T A  *  \  »UMMQMr  .  c  MN(jL  TD. RHCLUNr * » 

CAli.  PLT'  YTl?.40«r.»f..  •STTFfNPRr»  OJ  $T<?I  RUYIUN  N  .  I  %  *  »  .  0?  .  0  ,  > 

CALL  '•»L\U«R  f  4  ,?R  .*.*«>•  TK.  ?HT?  ..0  7.0.  » 

CALL  PL  Tr  XT<  ?  .4  .f  .f.  •^TRrNGTH  P  I'.  T  R  I  *»  UT  TUN  N.*%*  *  .07*0.  ) 

C  AL  L  ;’L.NIiM?.(  4  .?  .r,  .PO.  IS  •  ?Hll .  .07.0.  » 

CALL  PLNll.'MH  C  ».^»c  .F.FN  .NNCP.7HT4.0, 14.0.1 
CALL  UpL"7  (N1  ,NpLUK?.D?.n3.N4  1 
PPJNT  1  ON  NqL  OK  .NHL  OK 
C,f!  T(i  1 

f'  C«LL  1. 1  Nr  (HI  NPL  .YOATA.F  Itl’CC  •  YC  A  7  A  .N  .  4HNUNF  *  NHSUL  I CV  •  R'-*CLPSF  n  ) 
rp  Tm  1 

C  rilRMAT  NT  A  T  ph  rN  T  N  . . . . . . . 

101  r  HP  MAI  ( I  4  »?r  R  .  1 .  *T  1  .c  7. 7. 14  ) 

10?  n>kyAT  <1W]> 

JOY  FORMAT  (  1  *0.“4HnT,'PLa~r«-;rMT  LOftOTM'.  PFNULlR.^NX.Tta) 

104  rpRMAT  1 1  HO  .  I  4  ,  J  «  y  ,  1  1  HPF  AK  F(iPCrr.F7,0.  10MPFAK  H T :  . r 0 . 4 . 

1  OHM  A  v  .1?MYIFL0  FDPC.F  S  *  F  7. 0  .  *  l  M  V  1^1.  *>  PTS°=.r<,.4) 

1  O*  r  up  M#  T  (mC.MMRLtCK  AOQQFSS  .  IX  .  T \l>  »  *X .  ’HTt)  .  T1  f  1 
100  CALL  F N 0 R L T 
F  WO 


ppiNT  lO*.rK,vK  67 

PRINT  1  0°  •  f>»  VA  f>X 

1^1  riWMTHM  .  If  (FH*  **•*•) /?tw  5T  TFFMf  r>  v  (3  IS  TR  TPUT  TUN  I 
1  m  FORMA  T(  lHn’KX,  7H  U\  r  ♦="  n  P  M  1 
10^  r  f;RMA  T{  l*U»?4Xt9^UMM(jnH  ) 

1  CU  F  fit?  '■*  A  r  (  1 1!»  t  ?  H  X  »  I  >  MQ  T  G  MT  5KT  Wrn  I  1 

!  nr.  PflR  MA  7  (  1  M  ♦  «  ?4  X  »  1.1  ML  FF  T  SKPWfni 
1  OP  F nR  MA  T  (  IN*  «  ’(I  *  t  7-)  *7  MUOAL  I 

107  FORMAT  I  <IT3*  OyTItJW  H&HC,r  FROM » 1 PP I F. A «  JH  TU«1PFIK«»* 

iarxtiprjf  ,s> 


10F  FPRMA  T(  |  "MOni^TRT  PUTTl.N  *5  **• !*?!*•».  5*. 8MVAR TJ NCF . 1 PF 1 6  .  •» 

'  p  FT  UR  NJ  ™7<:  »<?AN*  lPf  15.F  »5X.  flHVARTANqr  ,  JPf  JF.p  » 


r  no 
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c usruut  jtir  rr  mk  (  m»  jk  ,  vk  » STirr  ) 

H  TMEN  SI  UN  STirrjNi 
CALL  A  Nil  NO  »NG1 
FKr POOOOO 0. /N 
cllMX~0. 
suMxrrn. 

so  T IJ  (1  ,  X  ,*  ,  7  #q  )  *  TK 
1  P-FK+S8RT  (  7,  *  VK  ) 
f-?.*EK-n 

no  ?  I r 1  .  m 

FALL  UN  r6-  PM  (A.H,XI 
^IIMXr  f,UMX*X 

su^x? :sum vr* i x*  x  I 
?  «  TIFF  U  1:  y 
on  Hi  n 

1  P "F  K*  SORT  (  F  .  *  VK  ) 

A “? .*FK  -p 

Of)  «  T  r  1  •  N 

CALL  IJM  IM  i  iO  (  A  «  P  »  X  | 

F|JNX'  SU^X  ♦  V 
5.UM  X?  rSUM  »’♦  I  X*  X  1 
•i  «TIFF(I  I r  X 
on  to  ii 

c  R=F  K* SORT  | ? . • VK I 
A"3.*FK-’.*R 
OP  f,  I;  i .  \ 

CALL  RSKF  V  t  A  «B  »  X  1 
<IJMX7SUWX»X 
F|JMX">rSlJMX?*(  X*  X  1 
f>  STIFC  U  |:  X 

on  tu  ii 

7  3-FK-»S0PT  (B.«VK  I 
ArX.«EK/*».-R/7. 

on  <?  1=1, M 
CALL  LSKF  v  (  A  *ft  t  X  ) 

SUMXr  SUMV*X 
sumx-=«su<x;>m  X*X> 

«  OTI  Ff  (I  ,T  x 

C'(?  T|.  11 

R  0  r  F  K ♦  S3  R  T  (  7  .  *  VK  ) 

A=?  .*CK-R 

no  If;  I  r]  .N 
CALL  'll  Mil  (  £  *P  t  X  I 
S  1IM  x :  SUMX  *■  X 
S(JMX^rSU’1X?*l  X*  X ) 

10  <tifF  IJ  i-  y 

11  DPINT  101 
XPA  R-Sljvy  /iv 

VAR  Xr  {F(JMX?/n  »-|  X-’AP.XR  AP  » 
on  TP  (1 7 . 1 T, m ,i c#  , f  IK 

ln  frimt  lo-» 

CP  Til  17 
1*  POT  NT  JOT 
on  TO  17 
IP  PRINT  10U 
CP  TU  17 
IS  POINT  lot 
or  Tp  17 
IF  o PI  NT  1  Of 
17  PRINT  107. A* ", NO 


"llRPilUTItir  r.fNS  INfjStVS.STRNTHtA  » 
n  r ^ r m  -i  tj»s»  s  tr  nt  ni  n> 

r  I/MX:  1. 
r  UM  X?  Ti. 

"r fl  Tff  H  .7t«c»7,q»,  re; 
i  n ?a  *r-  rjp  ti  i  ■>.* vs ) 

ESr  |  ft  *U  )/ ■>. 

nu  ?  I;  1*  N 

CALL  UM  IP  PM  I  A  «B  #  X  ) 

SUMX=SUMXl-X 
SUMX->TSijMX?*|  X*X> 

?  S  T»NT  HI  I)  rx 

r»r?  t  » t  n 

7  "TA  ♦VTPTI  ?<».♦  VS  I 
rrz  If  >/?. 
on  4  r:|iM 
CALL  UN  T  M  tlD  I  ft  » R  *  X  ) 

SUM  xr  SI* MX  »x 
r  (JM  x /’  rSl/MX?*|  X*X) 
u  S  TP  NT  M|  1 1  rx 
GO  Tn  n 

S  R  r*A  ♦<  OR  T  I  T  A  .  *  VS  I 
f  Sr  t?  ,*  !♦  A  1  /t  m 

nn  S  I:l.f| 

CALL  ^’SKr,,,  (A,9#X) 

* UM  X: SUMX  * X 
c  UM  X? tSUMX?4| X*  X  I 
*'  stpntuid-x 
Till  Tl?  11 

7  f>  :  A  +s  PR  Tl  1  fi ,  *  VS  I 
r  S:  (?.*a+o,)/*». 
nn  <?  r=i. m 

C«LL  LSKr  W  I  AtB  «  X  ) 

SUM  Xr 'UMX  *x 
SUMX?  TSUM  X?+ | X«  X  1 
S  «  TRNTMI  T)  :X 
SO  Tn  11 

9  R  r  A  ♦  «*.  OR  T  I  P,  *  VS  > 

T  S:  I A  *R  »/  •». 

no  hi  m  ,n 

CALL  HI  MM  o  1A.R.XI 
SUM  X:  S t J M X  ♦  X 
SUM  X?  r  Sl/M  X?*  1  X*  X  I 
1C  s  TRNTMI  II  rx 
n  X  Rfl  7r  SU  MX  /N 

VAPX: 1SUV X?/N 1- 1 XPAR^xn 4P ; 

ORINT  101 

cn  T !'  1  1  “  »  1  7  •  14  *  1  s  «  1  f  I  *  IS 
1  7  P°I  NT  1  0’’ 

C-U  T(i  l  ? 

1  7  PRJ  NT  1  O’ 
rn  Ti:  !  y 
!*»  PRINT  10*4 

nn  Ti*  17 

is  PRINT  1  OS 
rn  Ti.  i? 

If  P^IMT  in: 

17  PRI  NT  lfiM.a 

PRINT  io”.rs.vs 

P°IM7  lp.*»  ,  X>U;},VM-'X 


I 

w 

I 

l 

l 

l 


1  ?>!  F-tift  MA  t  (  ,*  ■*H0C,  T  *r  Nr  TH  ft  r  S  TRTRUTTUNJ 
1  0  1  f  I'R  Mf  T  (  1  M  ♦  •  ?u  *  #  7m  tjfg  ff  ) 

1  0  «  F  PR  MA  f  (  J  h  ♦  *  ?t»  X  *  «M  UN  I  M’.'D  m  ) 


1-0*1  F  (In  MS  T  (  1  H  +  ,  ->4X  ,  J 


hi?  I r  mt  SK  r  w**  I 


l  rtr  rnRM/i  T{  mlfft  skfwft,) 

1  r»r  r OR  MA  T  (  1  M  ♦  »  X  •  7M  >"»T  MUR  Aj.  ) 

107  P  np  ma  T  (  ■*«  HOP  T  ST  PI  QUTT  liN  R4N-r  r^uM.lPFlfc.^.XH  TCUIPFir 

rttuon  1  Hl"  T  T>  s  !>LE  MEAN,lprl£'«p*,>X*«HVAPlANrr,lPric.fl) 
r  np 


Reproduced  from 
best  available  copy. 


*.I|P  POUT  Ii^r  UNIFPMJA.R, 
P  -P  AN  UN  (  ■>  ) 

X  "  A  ♦('(-/»)  *  p 

-  r  r  u'>  n 

c  NO 


r.u«  Pl'UT  IN  -  UN  IMun  (  A  ,R,  X  > 

P  J  r  Rft  NUN(  P) 

P^Rfl  \'UN(  ?) 

*  *  C  . F  •  (  -1  ♦  ',+  ('5-A)*(R1*R7-i.)| 

P  T  T  (JP  N 

f  MD 


^U15  PlIUT  IN  f  ?SKF*I  5,0,x) 


p  1  -  RA  NUN( ?) 

PPrRA  NUN<  ,?> 

IFCP*.  3  F.P?)  r,(j  t ,)  J 

R-(  ?.  *C7,-»j  . 

r.l)  T  { *  7 

*P  1  )  -1  . 

X  *  0  .  r.  •  (  A  4  0  ♦  (  R  -  A  )  •  0  J 

P  r  T  Up  N 
r  ND 


^I'TRCUT  INf  l  F  KF  W  (  5,0.x) 
P 1 • Rft NUN( 7) 

NUN(  '?) 


rr(  PI  .Of.  •??)  r>U  T i)  J 
P-f  *P  1  )  -1. 
f?n  to  -> 


l  p-i?.*p?)-i. 

r  v  “0  *  (  .5,  0+  I  0 -a  )  •  P  ) 

P  F  T  UR  N 
F  ND 


*  UP  ROUT  INF  0  I  HUH  (  A  ,  0  ,  X  ) 
P I A PA NUM( ») 

P'T  Rfi  N[.'N(  ?) 

P"P  !♦».'- 1  . 

Tr(  p.  or  .T  .  )  r,,  T(,  2 

pa-  (j  ,  *  p  > 

r-n  Tc  ** 

l  p-i p 

*  “  0  . c  ♦  (  A  ♦  t  ♦  (  n  -  A  >  *  o ) 

P  F  T  Up  N 

run 
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APPENDIX  B 

TWO-DIMENSIONAL  PROGRAM 

The  following  pages  contain  a  listing  of  the  computer  program  for  the 
two-dimensional  model  including  four  required  subroutines,  the  data  deck 
for  the  Brazilian  test  mesh,  Figure  3.  3.  1,  and  the  plotting  program  with  its 
required  subroutine. 

The  plotting  routine  is  self-explanatory  with  the  documentation  in¬ 
cluded  in  the  subroutine.  The  subroutine  EZPLOT  was  written  for  the  IBM 
360/44  with  Calcomp  Plotter  and  probably  is  not  transferable  to  another 
machine.  It  is  Included  for  completeness  only. 

A  list  of  input  variables  (some  of  which  are  set  in  the  program  rather 
than  read  from  data  cards)  is  given  in  Table  B-l. 

The  program  takes  about  26  minutes  of  computing  time  on  an  IBM 
360/44  to  run  completely  through. 


1.  ISER 

2.  NP0IN 

3.  NELEM 

4.  NP 

5.  NIT 

6.  THICK 

7.  RAT 
B.  Xfl,  1) 

X(I,  2) 


TABLE  B-l 
INPUT  VARIABLES 

Serial  number  for  each 
Number  of  nodes 
Number  of  elements 

Code  number.  0  if  plane  strain.  Anything  else  if  plane 
stress. 

Number  of  nodes  to  which  an  external  load  is  applied. 
Thickness 

Ratio  of  compressive  to  tensile  strength 
X  coordinate  of  I  node 
Y  coordinate  of  1^  node 


73 


9*  U  (I)  ~  Initial  load  in  X  direction  on  l**1  node 

V(I)  -  Initial  load  in  Y  direction  on  I**1  node 

(Note  that  load  can  be  force  or  displacement  depending 
on  value  of  CGtoE) 

10.  CGtoE  (I)  -  0  if  U  and  V  are  forces,  1  if  V  is  a  force,  2  if  U  is  a  force, 

and  3  if  U  and  V  are  displacements.  If  code  is  greater 

than  three,  that  node  is  a  loaded  node  and  four  is  sub¬ 
tracted  from  its  code  before  the  code  is  used  by  the 
program. 


11. 

EPd,  1)  - 

Modulus  of  elasticity 

EP(1, 2)  - 

Poisson's  ratio 

EP(1. 3)  - 

A,]. 

Tensile  strength  for  I  element 

12. 

A  and  B 

Minimum  and  maximum  value,  respectively,  for  modulus 
of  elasticity,  defined  in  the  subroutine  CARDS 

13. 

C  and  D  - 

Minimum  and  maximum  value,  respectively,  for  tensile 
strength,  defined  in  the  subroutine  CARDS 
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